Propagation of a Huygens front through turbulent medium 
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The dynamics of a thin Huygens front propagating through turbulent medium is considered. 
A rigorous asymptotic expression for the effective velocity vf proportional to the front area is 
derived. The small-scale fluctuations of the front position are shown to be strongly intermittent. 
This intermittency plays a crucial role in establishing a steady state magnitude of the front velocity. 
The results are compared with experimental data. 
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The problem of propagation of a thin passive front 
(of flame, phase-transition, etc.) through turbulent flow 
has attracted a lot of attention since early forties [QJ|] 
when it was realized that velocity fluctuations tend to 
generate strongly convoluted front, thus dramatically in- 
creasing their area. In premixed combustion processes 
the flame front area is directly related to the speed of 
the front propagation. The first expression for the flame 
front velocity vp oc u rms , valid in the limit u rms 3> uq, 
where u rms is the root-mean-square velocity of the in- 
tegral eddy and uq is the laminar flame speed depend- 
ing on the details of chemical kinetics, was proposed by 
Schelkin Q (see for modern reviews on the theory 
of turbulent combustion) . The problem is also important 
for description of light propagation in the media with 
fluctuating dielectric constant, shock wave fronts etc. ||. 
The renormalization group approach 0, which included 
some theoretically unjustified steps, yielded the expres- 
sion, differing from the Schelkin result by a logarithmic 
factor, which agreed with experimental data [|s| — 1]~0|] in 
a wide range of parameter variation. Still, despite sub- 
stantial activity, no rigorous derivation of the front speed 
appeared and the question of the dependence of vp on 
both uq and the Reynolds number Re remained open. It 
is shown in this work that recent advances in the theory 
of a passive scalar, advected by turbulence, enable one to 
accurately account for the small-scale intermittency of a 
scalar field, crucial for description of the front fluctua- 
tions and derivation of the effective velocity vf. 

We consider a problem of propagation of a passive front 
though turbulent flow. The front can be described by 
the equation for a passive scalar PJTT) (so called "G"- 
equation) 



d t G + (vd T )G = u \d r Gl 



(1) 



where G(t,r) is a scalar field whose level surface, say 
G = 0, represents the thin front position. Statistics of 
turbulent velocity v is supposed to be known. The equa- 
tion (Q) with v = describes a front propagating with 
the constant speed uq (laminar front speed) normally to 
the local orientation of the front. For example, if the 
front at time t = 0, is defined at the x — y plane, it will 
propagate with a constant speed uq and constant area 



So = const in the z- direction. The role of the random 
field v is in generation of a strongly convoluted ( " wrin- 
kled") front with a substantially increased area St > So. 
When v — the mass of the reagents (fuel) consumed 
per unit time is: dm/dt = uqSq = const does not change 
in time. In general, dm/dt = uqSt = vpSo, where St 
is the area of the wrinkled front. This gives a definition 
of a turbulent or effective front velocity, vp = uqSt/So- 
Assuming that a steady (both vp and St do not depend 
on time after a long evolution) regime is realized, we in- 
troduce a new variable, G(t, r) = vp t — z + hit, r). Then, 
(111) reads: 



d t h + (vd r )h = v z + u y(d r h) + 1 - 2d z h - v F . (2) 

In the moving frame the fluctuations of the front position 
h are assumed to be in a statistically steady state, that 
fixes the value of the front speed vp, which in the turbu- 
lent regime discussed ( v F ^> u , when \d r h\ 2 ^> \d z h\ » 
1) is given by 



v F = u (\d r h\) - u (\Sh(r )\) /r . 



(3) 



Here, |<Wi(ro)| is a magnitude of the velocity difference 
at the scale ro where the "chemical" (uo-dependent) and 
advective contributions to (||) balance each other. Aver- 
aging over the turbulent velocity is assumed in (||). Since 
Urms 3> uo, ro <C i, where L is the scale of the turbu- 
lence source. This means that the scalar, injected at the 
scale L, is dissipated at the propagating front as a result 
of generation of very sharp cusps of the radius ro << L. 
Formation of such cusps was observed in numerical sim- 
ulations jl2|,[l3]]. Derivation of characteristic width r , 
magnitude |<5/i(ro)| of these cusps and, as a result, vp, is 
the goal of the theory. 

The presence of the two characteristic scales n and ro 
defines three possible flame regimes: 



A: L> r\ >ro] |B: L>r >li |C: L> 7? >rpl (4) 



Generation of the small-scale scalar fluctuations (direct 
cascade) 
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takes place in both inertial-convective, 
L > r > 77 and the dissipative-convective, rj > r > ro 
(which is valid in the A, C cases but not B) intervals. 
Thus, the problem is naturally divided into two: First, 
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we need to describe scalar (height of the flame brush) 
correlations in both convective ranges L > r ro- Once 
the solutions in the convective intervals are found, we will 
be ready to resolve the second and principal part of the 
problem: to calculate the value of the dissipative scale ro 
and, matching dissipative and convective intervals, find 
turbulent speed of the front, vp. Therefore, naturally, 
we are starting from the first task, considering all the 
regimes ( A — C) one after another. 

A. Batchelor (a pure viscous - convective) regime. The 
case A corresponds to a well- studied situation, first dis- 
cussed by Batchelor Q and developed further in [jl?] fl9| |. 
Without loss of generality and following we will con- 
sider the velocity difference to be Gaussian in the case 



\Sh r 



LyJJJfXyfln [L/r ] 



(9) 



D 



(5v?(t)6v?(t)) = - [25 a 



/3 r 2 



] exp 



(5) 



where r is turn-over time of the integral ( L-size) eddy. 
The pair correlation function of the scalar obeys the fa- 
mous logarithmic law in the convective interval of scales, 
L » r » r , ||§ 



{hih 2 ) = DL 2 In [L/r 12 ] / A, 



(6) 



where DL 2 in the nominator of the logarithmic prefactor 
describes the fluctuation of the " source" function v z while 
A stands for the Lyapunov exponent corresponding to the 
rate of Lagrangian stretching. Correlation between the 
source and convective terms in @ does not contribute to 
(0). Accounting for these correlations slightly modifies 
the higher-order moments generating subleading contri- 
butions and a mere renormalization of bare coefficients. 
The Lyapunov exponent A as a function of r, D was found 
for the two-dimensional version of the model (||) in [jl9| . 
Asymptotic of the large and small r were described ex- 
plicitly in p9[. Generally, the problem of finding A was 
reduced in |19| ] to a well-defined auxiliary quantum me- 
chanics which was easy to solve numerically. An interpo- 
lation formula for A, fitting well all the known asymptotic 
is: 



A = \/D/t tanh 



(7) 



This formula, valid in the space of arbitrary dimension- 
ality d > 2, holds up to 0{d) corrections. Statistics of 
the scalar fluctuations in the convective interval is shown 
to be Gaussian p9[| . Therefore, the typical fluctuation of 
the height h at r , estimated by the second moment (0) 



\5h r 



Ly/D/X. 



(8) 



The expression (J^) is derived without any spatial averag- 
ing over the large- scale ( ~ L) structures, always present 
in a real flow where all macroscopic characteristics, in- 
cluding the Lyapunov exponents A, are slightly modu- 
lated on the integral scale. Accounting for this spatial 
variation gives an estimate following directly from (^|): 



The modified regime A, accounting for the large- scale 
averaging, will be denoted hereafter by A'. 

B. Pure incrtial-convective range, fnertial-convective 
range is realized at the scales r S> r\ S> ro where one can- 
not neglect small-scale advection contribution . In the 
case of a general non-smooth velocity one finds a strongly 
intermittent behavior, manifested in the anomalous scal- 
ing of the scalar structure functions: 



S 2n (r) = ([h{ ri ) - h(T 2 )] 



2n 



^ r 12 i 



(10) 



valid at the separations r\ 2 <C L. The fundamental origin 
of the anomalous scaling, Q 2n < n(^ 2 , was discovered re- 
cently fHHH. it was understood that the anomalous ex- 
ponents originate from zero modes of the eddy-diffusivity 
operator: C, 2n are universal numbers, solely defined by 
the velocity statistics and independent on the proper- 
ties of the pumping term. The exponents were analyti- 
cally calculated for the case of the velocity field, rapidly 
varying in time, introduced by Kraichnan ]2(i[ in 
2 — ( 2 - ) and ^-expansions J||], respectively. 

An instanton approach (24| yields yet another large n ( 
n ^> d) asymptotic for Q 2n — > Coo(C2,^) when n — > oo. 
The constant Coo(C2)^) was explicitly calculated. The 
saturation of exponents C, n was predicted also in PH). 
The instanton consideration, applied to a general passive 
scalar problem, always results in the collapse of expo- 
nents (2n ~~ * Coo < oo with the asymptotic value Coo to 
be a complicated functional of the velocity field statis- 
tics. It can be easily understood: the 2n-th moment of 
the scalar can be represented as a path integral over 2n 
fluid particles. The dominant contribution into the high ( 
2n-th) order structure function originates from the most 
probable 2n- particle trajectory. In the incompressible 
world it is impossible to avoid a divergence of particles 
at least in one of the directions and it has been shown 



in 1 24 1 that contribution to the path integral from a sin- 
gle diverging trajectory is sufficient to cause saturation 
of the exponents Cn- This effect has a very strong and 
important influence on the parametric dependence of an 
effective front speed ut calculated below. 

C. Consecutive inertial-convective and dissipative- 
convective ranges. The description of the scalar corre- 
lations in the inertial-convective regime (ICR) does not 
deviate from the one considered above for the case B. 
However, in the dissipative-convective interval (DCR), 
L r\ r 3> ro, the behavior of the scalar is very 
different from that observed in the low-Reynolds- num- 
ber Batchelor regime when L sa r\. The crucial difference 
stems from the essential non-Gaussianity and intermit- 
tency of the scalar field on the velocity dissipation scale 
r\ where the solutions in both intervals have to match. 
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In this case the scale r\ is an integral scale for the DCR 
tq <C T <C r\ where powerful injection of all higher- order 
integrals of motion ( J dr h n ) takes place. The point is 
that a In — 2 moment defines a pumping for the next 2n- 
th order one. Thus, considering a generalization of the 
Kraichnan model for the case C( and neglecting the cor- 
relations between the source and convective terms, which 
can slightly renormalize some constants) one obtains the 
exact equations for the multi-point correlation functions 
^2n(ri,- • - ,r 2 „) = (h(ri)h(r 2 )...h(r 2n )) , 



LF 2n = V n (r 12 )F 2n ^ 2 (r 3 , • ■ • ,r 2n ) + perm. 



(11) 



where L 



-/C a/3 (r 



, V af) {r) = K af \L) - 
K a/3 (r) and K(r) ~ r 2 " 7 L» at, L > r > r], while 
K(r) r~j r t] r D at rj > r. Operator L has the —7 dimen- 
sionality (it scales as r -7 ) in the inertial-convective range 
while it is O(r ) in the dissipative-convective range. As 
a result, the 2n-th moment F 2n , explicitly depending on 
all 2n-th points, is dominated by the forced "logarith- 
mic" solution of (^Tj) in the DCR. In the ICR, however, 
it is a zero mode of the operator L which dominates the 
solution. Therefore, the 2n-th moment of the scalar dif- 
ference, S 2n (r) — (Sh 2n ) is estimated by 



S 2n {r) = a n L 2 In 



S 2n - 



(12) 



where a n are dimensionless n-dependent constants. On 
the another hand, at the scale 77 (|lj) this expression 
should match the anomalous structure functions (|iT 
from the upper (inertial-convective) interval giving: 



In" [r/r ]] 
In" [77/ r ] 



(13) 



Let us proceed now with the second task and estimate 
the value of the dissipative scale, 7*0, required for calcu- 
lation of the front velocity ut- The point of a crucial 



importance is a necessity to distinguish between the dif- 
ferent moments of the front fluctuations at the dissipative 
scale, Sh ro . Usually one estimates the typical fluctua- 
tions as a root-means-square value of the corresponding 
variable. This is fine in case of "normal scaling" or if 
the small-scale intermittency is not too strong (cases A 
and A'). However in the cases B and C intermittency 
in the convective intervals is extremely strong: estimates 
for a typical front fluctuation based on various moments, 



5h 



(n) 



(6h$ 



l/n 



are very different. As a result we get 



from (pU9| JlO| , |l3| ) the following n-dependent estimate at 
the dissipative scale rg 



\shi: 



L 



f L^D/X, A, 

I y/D/Xy/hl[L/r ], A', 

[ro/L}^\ B, (W) 



for the cases A,A',B and C respectively. On the other 
hand, the value of \Sh ro \ is defined by equilibration of the 
convective and "dissipative" terms in (0) giving: 



\5h%>\ 



v F T, T 



A- 1 In [L/r ], A, A' , 
L/v rms , B,C. 



(15) 



It shows that the linear dimension of the flame brush, 
calculated on the dissipative scale, is defined by the tur- 
bulent flame speed and the overall time of Lagrangian 
evolution which is a typical time for Lagrangian sepa- 
ration, initially equal to rg, to reach the integral scale. 
In the B, C regimes the overall time is tq- independent 
(note that the time of evolution from to r\ is neglected 
in comparison with one describing the inertial-convective 
stage of evolution, from n to L). Which number n is es- 
sential to define vp via ( |15|) is yet to be discussed. One 



should add (g) to ( 
one table we get 



Collecting all the relations in 



.4 

Ly/DX 



A' 

lVdx 



B 



F \n[L\/v ] ^ln[LA/u ] 

r fln[LX/v ] f\n[LX/v ] 



L 



fin 



C 

In" 1 / 2 



V. ( Vrms V 
L \ v J 



(16) 



where, /3 n = [1 — Ci + (n/n]~ ■ As was cited before, 
(2n/[2n\ — > and respectively, /3 n —> 1/[1 — Ci] > 0, 
at n — > 00. Since, intermittency of the scalar is grow- 
ing downscales from the scale of the pumping, the higher 
ratio of the integral scale to the dissipative one (of the 
turbulent velocity to the chemical one) the higher the 
moment-number n entering the actual stationary veloc- 
ity of the flame (which is defining the actual width of the 
front). The highest value of 6h$ is reached at n — > 00. 
That is to say that the stationary speed of the front is de- 
fined by such optimal and rare (instanton) configurations 



which give the dominant contribution into the highest 
moments of the scalar field h. The optimal configura- 
tions define the maximal spreading of the flame brush 
and control (mainly) the dissipation and the stationary 
velocity of the brush finally. The necessity to account for 
the highest moments of the Sh rd has a simple Lagrangian 
explanation. The point of a crucial importance is a fail- 
ure of 2n- particle description of the 2n-order correlation 
function of the scalar at the scales smaller than the dissi- 
pative one. This fact can be easily illustrated by a numer- 
ical procedure attempting to describe the initially plane 
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front in terms of collection of n Lagrangian particles |2q | 
each moving with velocity «on + v(x, t) where n is a unit 
vector in the z-direction. The initial ( t — 0) distance 
between these particles is Zo ~ 1/n. Very soon after be- 
ginning of the simulation the sharp cusps started to form 
and one had to add more and more particles to preserve 
continuity of the front. The development of these very 
sharp gradients, where dissipation takes place, drives the 
necessary number of particle n — » oo. We would like to 
reiterate a very important and profound peculiarity of 
the situation: the expression (|l6|), taken at the largest n, 
corresponds to a very rare optimal (instanton) configu- 
ration, responsible for the highest moments of the scalar 
difference. This means that to describe the small-scale 
( r <C r ) dynamics of the front height h(t;r ), one has 
to evaluate contributions from these configurations. We 
postpone the more accurate investigation ( which requires 
some further development of the instanton technique of 
pi| accounting for the entire dynamics of the front sur- 
face) for a future publication. 

One of the principle results of this work, given by Jl6| ) 
at n — ► oo, is that turbulent front speed strongly de- 
pends on the ratio T = ro/i], which is a novel dimen- 
sionless parameter of the problem. Therefore in order to 
compare experimental data with theoretical predictions 
we must first estimate T and determine the front propa- 
gation regime. Recent experimental studies of a passive 
front propagating in turbulent media were performed in 
the flows generated by vibrating grids, capillary waves, 
Taylor- Coucttc flow and Hell-Show cells. All experiments 
(for the turbulent regime, v rms 3> vq) corresponded to 
the case C . For example, in the vibrating grid experi- 
ments pfl] the values of u rms w 1 cm/ sec, L rs 1 — 10 cm, 
Re « HP - 10 3 and u /u rms w 1CT 2 - 1(T 3 . (In the 
case of the Kolmogorov turbulence one has an estima- 
tion r\ w 10LRe~ 3 / 4 , where Re is the Reynolds number 
and the factor c ~ 10 agrees with available experimental 
data). For the experimental conditions of || one gets 
r w 10 -3 — 10~ 4 . The expression for the front veloc- 
ity, following from the Table (|ilf), evaluated at ro corre- 
sponding to n — ► oo is: 



VF 



yjlln WL] + ln[U] V ln P\ ' 



(17) 



where U = u rms luQ. The rhs of ( |i7| ) (similar formula 
was derived in jj, see also |^7| , where the dynamic 
renormalization group has been used for evaluation of 
the turbulent flame speed vf) agrees very well with ex- 
perimental data on vp in a variety of turbulent flows 
in a wide range of variation of the dimensionless turbu- 
lent intensity U < 20 — 500. This universality can be 
readily understood: in all the experimental situations 
the 0(ln[i]/ L])- contribution is not large and can be 



neglected in comparison with the 0(ln [?7])-term. In a 
typical case U « 20 — 100 the transitional it rms -based 
Reynolds number is: Re c w 10 3 — 10 4 , which is very 
high. This explains the relatively broad applicability of 
the large U asymptotic of (|l7|) and therefore challenge 
a higher Re experiment to test the general structure of 
©■ 
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